library(tidyverse)
## -- Attaching packages ------------------------------------------------- tidyverse 1.2.1 --
## v ggplot2 3.2.1 v purrr 0.3.3
## v tibble 2.1.3 v dplyr 0.8.3
## v tidyr 1.0.0 v stringr 1.4.0
## v readr 1.3.1 v forcats 0.4.0
## -- Conflicts ---------------------------------------------------- tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()
library(magrittr)
##
## Attaching package: 'magrittr'
## The following object is masked from 'package:purrr':
##
## set_names
## The following object is masked from 'package:tidyr':
##
## extract
drivers <- feather::read_feather("drivers.feather")
incidents <- feather::read_feather("incidents.feather")
non_motorists <- feather::read_feather("non_motorists.feather")
combo <- feather::read_feather("combo.feather")
incident_locations <- incidents %>%
filter(!not_in_county) %>%
select(longitude, latitude)
kmeans_list <- lapply(1:15,function(x)kmeans(incident_locations, centers = x, nstart = 25, iter.max = 100))
## Warning: Quick-TRANSfer stage steps exceeded maximum (= 2852800)
centroidss <- data.frame(
centroids = 1:15,
sum_of_squares = sapply(1:15,function(x)kmeans_list[[x]]$tot.withinss)
)
centroidss %>%
ggplot(aes(x = centroids, y = sum_of_squares)) +
geom_line()

centroidss %>% knitr::kable()
| 1 |
805.93049 |
| 2 |
304.84448 |
| 3 |
213.07526 |
| 4 |
172.06474 |
| 5 |
139.81327 |
| 6 |
117.46724 |
| 7 |
99.63681 |
| 8 |
87.33599 |
| 9 |
75.62231 |
| 10 |
66.97816 |
| 11 |
60.92719 |
| 12 |
55.57872 |
| 13 |
50.99439 |
| 14 |
47.48665 |
| 15 |
44.29750 |
incident_clusters <- incident_locations
for(x in 1:15){
incident_clusters$cluster <-
as.factor(paste0("cluster", kmeans_list[[x]]$cluster))
temp <-incident_clusters %>%
ggplot(aes(x = longitude, y = latitude, color = cluster)) +
geom_point()
temp %>% print()
}















contingency_table <- function(data, var1, var2){
data %>%
group_by(!!enquo(var1), !!enquo(var2)) %>%
tally() %>%
ggplot() +
geom_bin2d(aes(x = !!enquo(var1), fill = log(n), y = !!enquo(var2))) +
geom_text(aes(x = !!enquo(var1), label = n, y = !!enquo(var2))) +
scale_fill_gradient(low = "#ece7f2",
high = "#a6bddb") +
theme(
panel.grid = element_blank(),
legend.position = "none",
panel.background = element_rect(fill = "white")
)
}
temp <-chisq.test(as.factor(incidents$weather), as.factor(incidents$acrs_report_type))
## Warning in chisq.test(as.factor(incidents$weather),
## as.factor(incidents$acrs_report_type)): Chi-squared approximation may be
## incorrect
temp %>% broom::tidy() %>% knitr::kable()
| 171.3451 |
0 |
24 |
Pearson’s Chi-squared test |
p-value 3.249696710^{-24}
contingency_table(incidents, acrs_report_type, weather)

Helmets
helmet_data <- non_motorists %>%
filter(as.logical(pedestrian_type_bicyclist)) %>%
select(safety_equipment,injury_severity) %>%
mutate(safety_equipment = c("no helmet","helmet")[1+as.integer(str_detect(tolower(safety_equipment), "helmet"))])
temp <-chisq.test(as.factor(helmet_data$safety_equipment), as.factor(helmet_data$injury_severity))
## Warning in chisq.test(as.factor(helmet_data$safety_equipment),
## as.factor(helmet_data$injury_severity)): Chi-squared approximation may be
## incorrect
temp %>% broom::tidy() %>% knitr::kable()
| 7.365619 |
0.117783 |
4 |
Pearson’s Chi-squared test |
contingency_table(helmet_data, safety_equipment,injury_severity)

Speed Limits
ordered_injury <- drivers %>%
count(injury_severity) %>%
arrange(desc(n)) %>%
pull(injury_severity)
speed_limit_data <-drivers %>%
select(speed_limit,injury_severity) %>%
mutate(
speed_limit = factor(c("Under 35", "35 - 45", "50+")[1 + findInterval(speed_limit, c(31, 46))], levels = c("Under 35", "35 - 45", "50+")),
injury_severity = factor(injury_severity, ordered_injury)
)
temp <-chisq.test(as.factor(speed_limit_data$speed_limit), as.factor(speed_limit_data$injury_severity))
## Warning in chisq.test(as.factor(speed_limit_data$speed_limit),
## as.factor(speed_limit_data$injury_severity)): Chi-squared approximation may be
## incorrect
temp %>% broom::tidy() %>% knitr::kable()
| 995.2863 |
0 |
8 |
Pearson’s Chi-squared test |
p-value 1.554337810^{-209}
contingency_table(speed_limit_data, speed_limit, injury_severity)

## Loading required package: Rcpp
## Loading required package: rlang
##
## Attaching package: 'rlang'
## The following object is masked from 'package:magrittr':
##
## set_names
## The following objects are masked from 'package:purrr':
##
## %@%, as_function, flatten, flatten_chr, flatten_dbl, flatten_int,
## flatten_lgl, flatten_raw, invoke, list_along, modify, prepend,
## splice
## Disabling daily seasonality. Run prophet with daily.seasonality=TRUE to override this.
##
## SAMPLING FOR MODEL 'prophet' NOW (CHAIN 1).
## Chain 1:
## Chain 1: Gradient evaluation took 0.002 seconds
## Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 20 seconds.
## Chain 1: Adjust your expectations accordingly!
## Chain 1:
## Chain 1:
## Chain 1: WARNING: There aren't enough warmup iterations to fit the
## Chain 1: three stages of adaptation as currently configured.
## Chain 1: Reducing each adaptation stage to 15%/75%/10% of
## Chain 1: the given number of warmup iterations:
## Chain 1: init_buffer = 7
## Chain 1: adapt_window = 38
## Chain 1: term_buffer = 5
## Chain 1:
## Chain 1: Iteration: 1 / 100 [ 1%] (Warmup)
## Chain 1: Iteration: 10 / 100 [ 10%] (Warmup)
## Chain 1: Iteration: 20 / 100 [ 20%] (Warmup)
## Chain 1: Iteration: 30 / 100 [ 30%] (Warmup)
## Chain 1: Iteration: 40 / 100 [ 40%] (Warmup)
## Chain 1: Iteration: 50 / 100 [ 50%] (Warmup)
## Chain 1: Iteration: 51 / 100 [ 51%] (Sampling)
## Chain 1: Iteration: 60 / 100 [ 60%] (Sampling)
## Chain 1: Iteration: 70 / 100 [ 70%] (Sampling)
## Chain 1: Iteration: 80 / 100 [ 80%] (Sampling)
## Chain 1: Iteration: 90 / 100 [ 90%] (Sampling)
## Chain 1: Iteration: 100 / 100 [100%] (Sampling)
## Chain 1:
## Chain 1: Elapsed Time: 10.751 seconds (Warm-up)
## Chain 1: 9.755 seconds (Sampling)
## Chain 1: 20.506 seconds (Total)
## Chain 1:
##
## SAMPLING FOR MODEL 'prophet' NOW (CHAIN 2).
## Chain 2:
## Chain 2: Gradient evaluation took 0 seconds
## Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Chain 2: Adjust your expectations accordingly!
## Chain 2:
## Chain 2:
## Chain 2: WARNING: There aren't enough warmup iterations to fit the
## Chain 2: three stages of adaptation as currently configured.
## Chain 2: Reducing each adaptation stage to 15%/75%/10% of
## Chain 2: the given number of warmup iterations:
## Chain 2: init_buffer = 7
## Chain 2: adapt_window = 38
## Chain 2: term_buffer = 5
## Chain 2:
## Chain 2: Iteration: 1 / 100 [ 1%] (Warmup)
## Chain 2: Iteration: 10 / 100 [ 10%] (Warmup)
## Chain 2: Iteration: 20 / 100 [ 20%] (Warmup)
## Chain 2: Iteration: 30 / 100 [ 30%] (Warmup)
## Chain 2: Iteration: 40 / 100 [ 40%] (Warmup)
## Chain 2: Iteration: 50 / 100 [ 50%] (Warmup)
## Chain 2: Iteration: 51 / 100 [ 51%] (Sampling)
## Chain 2: Iteration: 60 / 100 [ 60%] (Sampling)
## Chain 2: Iteration: 70 / 100 [ 70%] (Sampling)
## Chain 2: Iteration: 80 / 100 [ 80%] (Sampling)
## Chain 2: Iteration: 90 / 100 [ 90%] (Sampling)
## Chain 2: Iteration: 100 / 100 [100%] (Sampling)
## Chain 2:
## Chain 2: Elapsed Time: 8.017 seconds (Warm-up)
## Chain 2: 8.846 seconds (Sampling)
## Chain 2: 16.863 seconds (Total)
## Chain 2:
##
## SAMPLING FOR MODEL 'prophet' NOW (CHAIN 3).
## Chain 3:
## Chain 3: Gradient evaluation took 0.001 seconds
## Chain 3: 1000 transitions using 10 leapfrog steps per transition would take 10 seconds.
## Chain 3: Adjust your expectations accordingly!
## Chain 3:
## Chain 3:
## Chain 3: WARNING: There aren't enough warmup iterations to fit the
## Chain 3: three stages of adaptation as currently configured.
## Chain 3: Reducing each adaptation stage to 15%/75%/10% of
## Chain 3: the given number of warmup iterations:
## Chain 3: init_buffer = 7
## Chain 3: adapt_window = 38
## Chain 3: term_buffer = 5
## Chain 3:
## Chain 3: Iteration: 1 / 100 [ 1%] (Warmup)
## Chain 3: Iteration: 10 / 100 [ 10%] (Warmup)
## Chain 3: Iteration: 20 / 100 [ 20%] (Warmup)
## Chain 3: Iteration: 30 / 100 [ 30%] (Warmup)
## Chain 3: Iteration: 40 / 100 [ 40%] (Warmup)
## Chain 3: Iteration: 50 / 100 [ 50%] (Warmup)
## Chain 3: Iteration: 51 / 100 [ 51%] (Sampling)
## Chain 3: Iteration: 60 / 100 [ 60%] (Sampling)
## Chain 3: Iteration: 70 / 100 [ 70%] (Sampling)
## Chain 3: Iteration: 80 / 100 [ 80%] (Sampling)
## Chain 3: Iteration: 90 / 100 [ 90%] (Sampling)
## Chain 3: Iteration: 100 / 100 [100%] (Sampling)
## Chain 3:
## Chain 3: Elapsed Time: 8.967 seconds (Warm-up)
## Chain 3: 10.215 seconds (Sampling)
## Chain 3: 19.182 seconds (Total)
## Chain 3:
##
## SAMPLING FOR MODEL 'prophet' NOW (CHAIN 4).
## Chain 4:
## Chain 4: Gradient evaluation took 0.001 seconds
## Chain 4: 1000 transitions using 10 leapfrog steps per transition would take 10 seconds.
## Chain 4: Adjust your expectations accordingly!
## Chain 4:
## Chain 4:
## Chain 4: WARNING: There aren't enough warmup iterations to fit the
## Chain 4: three stages of adaptation as currently configured.
## Chain 4: Reducing each adaptation stage to 15%/75%/10% of
## Chain 4: the given number of warmup iterations:
## Chain 4: init_buffer = 7
## Chain 4: adapt_window = 38
## Chain 4: term_buffer = 5
## Chain 4:
## Chain 4: Iteration: 1 / 100 [ 1%] (Warmup)
## Chain 4: Iteration: 10 / 100 [ 10%] (Warmup)
## Chain 4: Iteration: 20 / 100 [ 20%] (Warmup)
## Chain 4: Iteration: 30 / 100 [ 30%] (Warmup)
## Chain 4: Iteration: 40 / 100 [ 40%] (Warmup)
## Chain 4: Iteration: 50 / 100 [ 50%] (Warmup)
## Chain 4: Iteration: 51 / 100 [ 51%] (Sampling)
## Chain 4: Iteration: 60 / 100 [ 60%] (Sampling)
## Chain 4: Iteration: 70 / 100 [ 70%] (Sampling)
## Chain 4: Iteration: 80 / 100 [ 80%] (Sampling)
## Chain 4: Iteration: 90 / 100 [ 90%] (Sampling)
## Chain 4: Iteration: 100 / 100 [100%] (Sampling)
## Chain 4:
## Chain 4: Elapsed Time: 9.467 seconds (Warm-up)
## Chain 4: 9.49 seconds (Sampling)
## Chain 4: 18.957 seconds (Total)
## Chain 4:
## Warning: The largest R-hat is 1.11, indicating chains have not mixed.
## Running the chains for more iterations may help. See
## http://mc-stan.org/misc/warnings.html#r-hat
## Warning: Bulk Effective Samples Size (ESS) is too low, indicating posterior means and medians may be unreliable.
## Running the chains for more iterations may help. See
## http://mc-stan.org/misc/warnings.html#bulk-ess
## Warning: Tail Effective Samples Size (ESS) is too low, indicating posterior variances and tail quantiles may be unreliable.
## Running the chains for more iterations may help. See
## http://mc-stan.org/misc/warnings.html#tail-ess
plot(total_model, forecast)

prophet_plot_components(total_model, forecast)

## Disabling daily seasonality. Run prophet with daily.seasonality=TRUE to override this.
##
## SAMPLING FOR MODEL 'prophet' NOW (CHAIN 1).
## Chain 1:
## Chain 1: Gradient evaluation took 0 seconds
## Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Chain 1: Adjust your expectations accordingly!
## Chain 1:
## Chain 1:
## Chain 1: WARNING: There aren't enough warmup iterations to fit the
## Chain 1: three stages of adaptation as currently configured.
## Chain 1: Reducing each adaptation stage to 15%/75%/10% of
## Chain 1: the given number of warmup iterations:
## Chain 1: init_buffer = 7
## Chain 1: adapt_window = 38
## Chain 1: term_buffer = 5
## Chain 1:
## Chain 1: Iteration: 1 / 100 [ 1%] (Warmup)
## Chain 1: Iteration: 10 / 100 [ 10%] (Warmup)
## Chain 1: Iteration: 20 / 100 [ 20%] (Warmup)
## Chain 1: Iteration: 30 / 100 [ 30%] (Warmup)
## Chain 1: Iteration: 40 / 100 [ 40%] (Warmup)
## Chain 1: Iteration: 50 / 100 [ 50%] (Warmup)
## Chain 1: Iteration: 51 / 100 [ 51%] (Sampling)
## Chain 1: Iteration: 60 / 100 [ 60%] (Sampling)
## Chain 1: Iteration: 70 / 100 [ 70%] (Sampling)
## Chain 1: Iteration: 80 / 100 [ 80%] (Sampling)
## Chain 1: Iteration: 90 / 100 [ 90%] (Sampling)
## Chain 1: Iteration: 100 / 100 [100%] (Sampling)
## Chain 1:
## Chain 1: Elapsed Time: 5.801 seconds (Warm-up)
## Chain 1: 5.54 seconds (Sampling)
## Chain 1: 11.341 seconds (Total)
## Chain 1:
##
## SAMPLING FOR MODEL 'prophet' NOW (CHAIN 2).
## Chain 2:
## Chain 2: Gradient evaluation took 0 seconds
## Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Chain 2: Adjust your expectations accordingly!
## Chain 2:
## Chain 2:
## Chain 2: WARNING: There aren't enough warmup iterations to fit the
## Chain 2: three stages of adaptation as currently configured.
## Chain 2: Reducing each adaptation stage to 15%/75%/10% of
## Chain 2: the given number of warmup iterations:
## Chain 2: init_buffer = 7
## Chain 2: adapt_window = 38
## Chain 2: term_buffer = 5
## Chain 2:
## Chain 2: Iteration: 1 / 100 [ 1%] (Warmup)
## Chain 2: Iteration: 10 / 100 [ 10%] (Warmup)
## Chain 2: Iteration: 20 / 100 [ 20%] (Warmup)
## Chain 2: Iteration: 30 / 100 [ 30%] (Warmup)
## Chain 2: Iteration: 40 / 100 [ 40%] (Warmup)
## Chain 2: Iteration: 50 / 100 [ 50%] (Warmup)
## Chain 2: Iteration: 51 / 100 [ 51%] (Sampling)
## Chain 2: Iteration: 60 / 100 [ 60%] (Sampling)
## Chain 2: Iteration: 70 / 100 [ 70%] (Sampling)
## Chain 2: Iteration: 80 / 100 [ 80%] (Sampling)
## Chain 2: Iteration: 90 / 100 [ 90%] (Sampling)
## Chain 2: Iteration: 100 / 100 [100%] (Sampling)
## Chain 2:
## Chain 2: Elapsed Time: 6.701 seconds (Warm-up)
## Chain 2: 6.932 seconds (Sampling)
## Chain 2: 13.633 seconds (Total)
## Chain 2:
##
## SAMPLING FOR MODEL 'prophet' NOW (CHAIN 3).
## Chain 3:
## Chain 3: Gradient evaluation took 0 seconds
## Chain 3: 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Chain 3: Adjust your expectations accordingly!
## Chain 3:
## Chain 3:
## Chain 3: WARNING: There aren't enough warmup iterations to fit the
## Chain 3: three stages of adaptation as currently configured.
## Chain 3: Reducing each adaptation stage to 15%/75%/10% of
## Chain 3: the given number of warmup iterations:
## Chain 3: init_buffer = 7
## Chain 3: adapt_window = 38
## Chain 3: term_buffer = 5
## Chain 3:
## Chain 3: Iteration: 1 / 100 [ 1%] (Warmup)
## Chain 3: Iteration: 10 / 100 [ 10%] (Warmup)
## Chain 3: Iteration: 20 / 100 [ 20%] (Warmup)
## Chain 3: Iteration: 30 / 100 [ 30%] (Warmup)
## Chain 3: Iteration: 40 / 100 [ 40%] (Warmup)
## Chain 3: Iteration: 50 / 100 [ 50%] (Warmup)
## Chain 3: Iteration: 51 / 100 [ 51%] (Sampling)
## Chain 3: Iteration: 60 / 100 [ 60%] (Sampling)
## Chain 3: Iteration: 70 / 100 [ 70%] (Sampling)
## Chain 3: Iteration: 80 / 100 [ 80%] (Sampling)
## Chain 3: Iteration: 90 / 100 [ 90%] (Sampling)
## Chain 3: Iteration: 100 / 100 [100%] (Sampling)
## Chain 3:
## Chain 3: Elapsed Time: 6.286 seconds (Warm-up)
## Chain 3: 6.644 seconds (Sampling)
## Chain 3: 12.93 seconds (Total)
## Chain 3:
##
## SAMPLING FOR MODEL 'prophet' NOW (CHAIN 4).
## Chain 4:
## Chain 4: Gradient evaluation took 0 seconds
## Chain 4: 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
## Chain 4: Adjust your expectations accordingly!
## Chain 4:
## Chain 4:
## Chain 4: WARNING: There aren't enough warmup iterations to fit the
## Chain 4: three stages of adaptation as currently configured.
## Chain 4: Reducing each adaptation stage to 15%/75%/10% of
## Chain 4: the given number of warmup iterations:
## Chain 4: init_buffer = 7
## Chain 4: adapt_window = 38
## Chain 4: term_buffer = 5
## Chain 4:
## Chain 4: Iteration: 1 / 100 [ 1%] (Warmup)
## Chain 4: Iteration: 10 / 100 [ 10%] (Warmup)
## Chain 4: Iteration: 20 / 100 [ 20%] (Warmup)
## Chain 4: Iteration: 30 / 100 [ 30%] (Warmup)
## Chain 4: Iteration: 40 / 100 [ 40%] (Warmup)
## Chain 4: Iteration: 50 / 100 [ 50%] (Warmup)
## Chain 4: Iteration: 51 / 100 [ 51%] (Sampling)
## Chain 4: Iteration: 60 / 100 [ 60%] (Sampling)
## Chain 4: Iteration: 70 / 100 [ 70%] (Sampling)
## Chain 4: Iteration: 80 / 100 [ 80%] (Sampling)
## Chain 4: Iteration: 90 / 100 [ 90%] (Sampling)
## Chain 4: Iteration: 100 / 100 [100%] (Sampling)
## Chain 4:
## Chain 4: Elapsed Time: 6.918 seconds (Warm-up)
## Chain 4: 10.483 seconds (Sampling)
## Chain 4: 17.401 seconds (Total)
## Chain 4:
## Warning: There were 11 transitions after warmup that exceeded the maximum treedepth. Increase max_treedepth above 10. See
## http://mc-stan.org/misc/warnings.html#maximum-treedepth-exceeded
## Warning: Examine the pairs() plot to diagnose sampling problems
## Warning: The largest R-hat is 2.51, indicating chains have not mixed.
## Running the chains for more iterations may help. See
## http://mc-stan.org/misc/warnings.html#r-hat
## Warning: Bulk Effective Samples Size (ESS) is too low, indicating posterior means and medians may be unreliable.
## Running the chains for more iterations may help. See
## http://mc-stan.org/misc/warnings.html#bulk-ess
## Warning: Tail Effective Samples Size (ESS) is too low, indicating posterior variances and tail quantiles may be unreliable.
## Running the chains for more iterations may help. See
## http://mc-stan.org/misc/warnings.html#tail-ess
plot(alcohol_model, forecast)

prophet_plot_components(alcohol_model, forecast)
